######实践六：R与转录组分析
rna<-read.csv("mRNA_exprSet.csv")
library(DESeq2)
sample_info <- data.frame(  
  sampleID = colnames(rna)[-1],   #colnames(expr_data)[-1] 提取了rna数据框中除了第一列之外的所有列名。
  condition = ifelse(grepl("01$", colnames(rna)[-1]), "tumor", "normal"),  
  #使用正则表达式检查每个样本ID的后两位是否为 "01"。如果是，返回 TRUE,为“tumor'；否则，返回 FALSE,为”normal"。
  stringsAsFactors = TRUE)#用于指定是否应该将字符串自动转换为因子
#使用DESeqDataSetFromMatrix函数来构建一个DESeqDataSet对象。
dds <- DESeqDataSetFromMatrix (countData = as.matrix(rna[,-1]),  
                              colData = sample_info,  
                              design = ~ condition)
#使用DESeq函数来执行差异表达分析
dds2<- DESeq(dds)
results_table <- results(dds2)
results_table<-as.data.frame(results_table)
sig_genes <- data.frame(results_table$pvalue,results_table$log2FoldChange)
colnames(sig_genes)<-c("pvalue","logFC")
#绘制火山图
library("ggplot2")
FC1<-log(1.2)
FC2<--log(1.2)
up_regulated <- sig_genes$pvalue < 0.05 & sig_genes$logFC > FC1 
down_regulated <-  sig_genes$pvalue < 0.05 & sig_genes$logFC < FC2
sig_genes$Group<-ifelse(up_regulated,"up",ifelse(down_regulated,"down","Not Significant"))
ggplot(sig_genes, aes(x = logFC, y = -log10(pvalue), colour = Group))+
  geom_point(size =4, shape = 20, stroke = 0.5)+
  #控制最人气泡和最小气泡，调节气泡相对大小
  scale_size(limits = c(2,16))+
  ##设置颜色
  #scale_fill_manual(values = c("#fe0000","#13fc00","#bdbdbd"))+
  scale_color_manual(values=c('steelblue','gray','brown'))+
  ylab('-log10 (pvalue)')+
  xlab('log2 (FoldChange)')+
  ## 增加横竖线条
  geom_hline(yintercept=-log10(0.05), col="darkgray", linetype="dashed") +  
  geom_vline(xintercept=c(-log2(1.2), log2(1.2)), col="darkgray", linetype="dashed")

